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Solid state devices for quantum bit computation (qubits) are not perfect isolated two-level systems, 
since additional higher energy levels always exist. One example is the Josephson flux qubit, which 
consists on a mesoscopic SQUID loop with three Josephson junctions operated at or near a magnetic 
flux of half quantum. We study intrinsic leakage effects, i.e., direct transitions from the allowed 
qubit states to higher excited states of the system during the application of pulses for quantum 
computation operations. The system is started in the ground state and rf- magnetic field pulses are 
applied at the qubit resonant frequency with pulse intensity f p . A perturbative calculation of the 
average leakage for small f p is performed for this case, obtaining that the leakage is quadratic in 
f p , and that it depends mainly on the matrix elements of the supercurrent. Numerical simulations 
of the time dependent Schrodinger equation corresponding to the full Hamiltonian of this device 
were also performed. From the simulations we obtain the value of f p above which the two-level 
approximation breaks down, and we estimate the maximum Rabi frequency that can be achieved. 
We study the leakage as a function of the ratio a among the Josephson energies of the junctions of 
the device, obtaining the best value for minimum leakage (a ~ 0.85). The effects of flux noise on 
the leakage are also discussed. 

PACS numbers: 03.67.Lx,82.25.Cp,74.50.+r 



I. INTRODUCTION 

In the last decade different devices of mesoscopic 
Josephson junctions have been studied experimentally as 
candidates for the design of quantum bits (qubits) i - — 
A large effort has been devoted to succeed in the coher- 
ent manipulation of their quantum states in a control- 
lable wayji^^ One of the superconductig qubit devices 
that has been studied in the last years is the Josephson 
flux qubit,— which consists on a mesoscopic SQUID 
loop with three Josephson junctions operated at or near 
a magnetic flux of half quantum. 

Real qubit devices, however, are not perfect isolated 
two-level systems. First, coupling to the external en- 
vironment induces relaxation and dephasing.^ Second, 
additional higher energy levels always exist in solid 
state devices. Therefore leakage effects, i.e., transitions 
from the allowed qubit states to higher excited states 
of the system can occur during quantum computation 
operations-^ 2 Indirect leakage to the higher energy lev- 
els produced through the interaction with the environ- 
ment has been studied in some casesi 13 ! 14 Even neglecting 
the interaction with the external environment, intrinsic 
leakage can occur due to direct transitions outside the 
computational subspace during the application of pulses 
for computational operations! 11 i 15 i 16 Due to the impor- 
tance of minimizing the gate errors due to leakage, several 
optimization strategies to compensate leakage, based on 
varying the pulse shapes and pulse sequences, have been 
studied recently— ~— 

Furthermore, the study of the multilevel dynamics of 
qubit devices has become of interest by itself in the 
last years— r— The consideration of the superconducting 
qubit devices as artificial atoms, has lead to the study 



of the dynamic effects of their level structure beyond 
the lowest energy levels. Effects of strong drive ampli- 
tudes on Rabi oscillations have been studied—"— Driv- 
ing the flux qubit with large amplitude harmonic excita- 
tions have also revealed the higher energy levels through 
Landau-Zener-Stuckelberg transitions— — Mach-Zender 
interferometr y 25 ' 26 and amplitude spectroscopy^! have 
been the subject of recent studies of the flux qubit as 
an artificial atom. Moreover, it has also been pointed 
out that the high energy level structure of the Josephson 
flux qubit should show quantum signatures of classical 
chaos— 

For quantum computation applications one wants to 
maximize the number of quantum bit operations be- 
fore gate errors become important. In superconducting 
qubits, long pulses are limited by the decoherence due 
to the environment and short pulses by leakage out of 
the qubit computational subspace. In order to maximize 
the number of quantum bit operations, one has to max- 
imize the ratio ideph/iop, where idcph is the dephasing 
time and t op is the time scale for a single quantum op- 
eration. The main approach has been to improve the 
design of the qubit devices to increase their coherence 
time. In the case of the Josephson flux qubit there has 
been an important progress in increasing idcph, from the 
early experiments by Chiorescu et al& with idcph ~ 20ns 
to recent experiments that report tdoph ~ 0.5 — 2/is 4 Pro- 
vided that one has succeeded to achieve a ideph as large 
as possible for a given device, the following approach is to 
reduce t op . It is in this later case (t op <C ^deph) when the 
effect of leakage is relevant. The usual strategy for quan- 
tum operations is to drive the qubit with a periodic pulse 
of intensity f p at a resonant frequency hiu r = E\ — Eq, 
with Ei — Eq the energy difference between the two qubit 



states. In this case, the time scale for quantum comput- 
ing operations, t op , is proportional to the period of Rabi 
oscillations, t op ~ Tr. Since the time Tr depends on 
the pulse strength as Tr ~ to reduce i op , one has 

to increase f p . Along this line of reasoning, the following 
questions will be addressed here: (i) how much is possible 
to increase f p before leakage effects become important, 
and (ii) for which circuit parameters of the Josephson 
flux qubit the intrinsic leakage is minimum. 

To this end, in this work we will study the quantum dy- 
namics of the Josephson flux qubit solving its time depen- 
dent Schrodinger equation considering the full hamilto- 
nian of the system. Since we are interested in time scales 
such that i op <C idcptn the interaction with the environ- 
ment will be neglected, and we focus on the calculation 
of the amount of intrinsic leakage. We will study the case 
when the qubit is driven by an rf pulse in the magnetic 
field that is resonant with hio r = E\ — Eq . The amount of 
leakage as a function of the pulse strength f p will be cal- 
culated both perturbatively for small f p and numerically 
for arbitrary values of f p . The paper is organized as fol- 
lows. In Sec. II we introduce the model hamiltonian and 
equations for the Josephson flux qubit. In Sec. Ill we 
present a perturbative calculation of the leakage in the 
case of an harmonic resonant drive. In Sec. IV we present 
our numerical results for the time dependent Schrodinger 
equation, calculating the amount of leakage as a function 
of f p . In Sec.V we show results on the dependence of te 
leakage with different circuit parameters. In particuar, 
the optimum value of the circuit parameter a (the ratio 
among the Josephson couplings of the junctions of the 
qubit) will be computed. In Sec. VI we analyze the ef- 
fect of an small amount of noise in the results shown in 
the previous section. Finally, Sec. VII contains a sum- 
mary and a discussion of the most relevant points of our 
findings. 



II. MODEL FOR THE DEVICE FOR THE 
JOSEPHSON FLUX QUBIT 

The device used for the Josephson flux qubit 2 consists 
on a superconducting ring with three Josephson junctions 
enclosing a magnetic flux $ = f$o, with $o = h/2e, see 
FiglU 

The junctions have gauge invariant phase differences 
defined as (pi, ip 2 and ip 3 , respectively, with the sign con- 
vention corresponding to the directions indicated by the 
arrows in FigfTJ Typically the circuit inductance can 
be neglected and the phase difference of the third junc- 
tion is: <p3 = —ipi + ip2 — 27r/. Therefore the system 
can be described with two dynamical variables: <pi,(p2 ■ 
The circuits that are used for the Josephson flux qubit 
have two of the junctions with the same coupling energy, 
Ej t i = E Jj2 — Ej, and capacitance, C\ = C 2 = C, while 
the third junction has smaller coupling Ej 3 — aEj and 
capacitance C3 = aC, with 0.5 < a < 1. The above 
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FIG. 1: Circuit for the Device for the Josephson Flux Qubit 
as described in the text. Josephson junctions 1 and 2 have 
Josephson energy Ej and capacitance C, and junction 3 has 
Josephson energy and capacitance a times smaller. The ar- 
rows indicate the sign convention for defining the gauge in- 
variant phase differences. The circuit encloses a magnetic flux 
$ = /$ . 



considerations lead to the Hamiltonian 2 - 

H = ~P T M- 1 P+EjV(<p) ) (1) 

where the two-dimensional coordinate is <p = ((pi,<p2)- 
The potential energy term is given by the Josephson en- 
ergy of the circuit and, in units of Ej, is: 

V(<p) = 2+a—cos(pi — cos(p2—acos(2Trf+(pi~<p 2 ) . (2) 

The kinetic energy term is given by the electrostatic en- 
ergy of the circuit, where the two-dimensional momen- 
tum is 

P=(P 1 ,P 2 ) = M .^, 

and M is an effective mass tensor determined by the ca- 
pacitances of the circuit, 

with 



We neglected in M the on-site capacitances C g (typically 
Cg/C ~ 1CT 2 - l(T 3 < 1). The system modelled with 
Eqs. CE}-© * s analogous to a particle with anysotropic 
mass M in a two-dimensional periodic potential V(ip). 

In typical junctions, the Josephson energy scale, Ej, 
is much larger than the electrostatic energy of electrons, 
Ec = e 2 /2C, and the system is in a classical regime. On 
the other hand, mesoscopic junctions (with small area) 
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can have Ej ~ Eq, and quantum fluctuations become 
important. In this case, the quantum momentum opera- 
tor is defined as 



P -> P = -iW„ = 



After replacing the above defined operator P in the 
Hamiltonian of Eq. (JXJ) , the eigenvalue Schrodinger equa- 
tion becomes 



* v {0) = E^ v {0) , (3) 



where we normalized energy by Ej and momentum by 
h/ ^8Ec/Ej. We defined in Eq.Q the parameter r\ = 
J 8Ec /Ej which plays the role of an effective h~& Typ- 
ical flux qubit experiments have values of a in the range 
0.6 - 0.9 and r\ in the range 0.1 - 

In this work, we will study the quantum dynamics of 
the Josephson flux qubit. Therefore, we solve the time- 
dependent Schrodinger equation which, with the same 
normalization as above, is given by 



.gg(g) 
1 dt 



V{0) 



*(<?) , (4) 



where we normalized time by tj = h/Ej. 

We integrate numerically Eq. (Q} with a second or- 
der split-operator algorithm, 29 using a discretization grid 
of Atp = 2tt/I28 and At = O.ltj. We use 27r-periodic 
boundary conditions on <p = ((pi, (p^j. Eigenstates 
and eigenenergies Ei are also calculated by numerical di- 
agonalization of Eq. ([3]), with the same discretization grid 
and boundary conditions. In what follows we will con- 
sider the case of r\ = 0.48 (i.e., Ej/Ec — 35), which 
corresponds to the experiment of Chiorescu et al.. 6 



III. PERTURB ATIVE CALCULATION OF THE 
INTRINSIC LEAKAGE 

In quantum computation implementations 2 ' 6,7 the 
Josephson flux qubit is operated at magnetic fields near 
the half-flux quantum, / = 1/2 + Sf, with Sf <C I. 
For values of a > 1/2, the potential of Eq.@ has two 
well defined minima. At the optimal operation point 
fo = 1/2, the two lowest energy eigenstates (l^o) an d 
l^i)) are symmetric and antisymmetric superpositions 
of two states corresponding to macroscopic persistent 
currents of opposite sign. A two-level truncation of the 
Hilbert space is usually performed. 2 In the subspace ex- 
panded by |^o) an d l^i): the hamiltonian of Eq. (p} is 
reduced to 



= ~2 <Jz 2 <7x ' 



(5) 



where H is written in the qubit basis defined by |0) = 
(1*0) + |*i»A/2 and |1> = (|1<o) - |*i»/\/2. Here 



A = Ei — Eq is the two-level spliting at / = 1/2, which 
increases strongly with a, and e = AnaEjSoiSf (consid- 
ering that Sf < 1), with S Q1 = (* | sin(7r + ^i-^ 2 )|*i)- 
(For typical values of a and rj, one has Soi ~ 0.8). 
Most experiments control the system varying the mag- 
netic field Sf = f — 1/2. Recently it has been shown 
experimentally that is also possible to manipulate the 
value of A by controlling a, replacing the third junction 
by an additional SQUID loop^iS 

Most of the experiments on the flux qubit study the 
possibility of single bit quantum operations by driving 
the qubit with a resonant pulse in the magnetic flux with 
Sf(t) = f p sm(uj r t), at the resonant frequency huj r = Ei — 
Eq. If the system is started in the ground state \^>o), 
after the pulse is applied during a time interval r the 
populations of the ground state and the excited state are 



Po = |(*(r)|* )| 2 
Pi = K*(r)|*!)| 2 



cos 2 (tt fl T/2) , 
sin 2 (n fl r/2) , 



with the Rabi frequency MIr — e p /2 and e p ~ 
AiTaEjSoifp. This result is usually obtained in the ro- 
tating wave approximation (see below). 

In order to check the goodness of the two-level ap- 
proximation we are going to evaluate, perturbatively, 
the population of the higher energy levels when a pulse 
8 fit) = f p sm(ujt) is applied, with f p <C 1. We calculate 
the leakage outside of the quantum computational space 
spanned by the two lowest levels as 

oo oo 



where the \^f n ) are the eigenstates at fo = 1/2. 

We now write the hamiltonian of Eq. (JT|) as T~L = 
n(fo) + W(5f(t)) with f = 1/2 and 

W(Sf(t)) = aEj sin[27r<5/(t)] sin(7r + fa - fa) + 

2aEj sin 2 [7rd f(t)} cos(7r + fa - fa) , (6) 

for small f p we have 

W{6f(t)) w 2nfpaE.j sin(wt) sin(7r + fa - fa) . 

For the perturbative calculations we use the fact that 
the first and second states interact strongly with each 
other but only weakly with higher states. In this approx- 
imation we solve (see for example Ref 30) 



ih 



dc Q (t) 
dt 



= W 00 {t)c (t) + Woi(t)e Wolt a(t) , (7) 



dt 



ih 



dc n (t) 
dt 



W n0 (ty^ ot co(t) + WmWe^ait) , (9) 
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with Wij = (^i\W\^j) and uj 13 = (E l - E 3 )/h. 
We rewrite equations ([7]) and © in the form 

= n oo sin (ujt)c (t) +rtsin(ujt)e- lullot c 1 (t) , 



d Cl {t) 



(10) 



rjsin(^)e^ 10 *c (i) + fin sin (wt)ci(i) , (11) 



where we defined ftfiy = 27r f p aSijEj with Sly = 
(^sh^TT + 0i - ^ 2 )|%), and fi 10 = fi i = fi (for Sy 
real). In order to solve these equations we make the fol- 
lowing change of variables^ 



b 3 (t) = Ci (t)e i(n « /w)cos(u,t) 



then 



= r!sin(^)e"' i( " 10t+iycosM)) 6i(i) 



and 



where i/ = (fin — fioo)/^- Using the relation 

oo 

ixcos9 \ y -k t / \ mf? 

e = 2^ z Jfc(a;)e 

we obtain 



(12) 



(13) 



h— — oo 



.S6o(t) fi . fc . , . 



dt 2i 



-i(u>io — (k+l)co)t 



i(o;io-(fc-l)w)tl i ( 14 ) 



.0&l(t) fi ^ . fc _ . . 



dt 



2i 



O i(bjio + (k—l)u>)t 



h(t) , (15) 



where the Jk are Bessel functions. When uj ~ wio we 
can use the rotating wave approximation (RWA) 31 ' 32 and 
negelect the highly off resonant terms obtaining: 



db Q (t) 
dt 



and 



(16) 



(17) 



where 7 = (f2/2)(Jo(i/) + J2{v)). In the exact resonance 
case (lu — u>io) the problem has analytic solution. 

db Q (t) db.it) 

—Qf~ = -lbi(t); = 760 (*) . (18) 

If the system was initially in the ground state, we obtain 



c {t) = e ^ooM(i-cosM)) cos ( 7i ) 



(19) 



and 



d(t) =e ^oo/u ) -(n 11 /u ) )cos( W t)) sin (' 7t ^ ^qj 

We see that |c (i)| 2 = cos 2 (7i) and | Cl (t)| 2 = sin 2 ( 7 i), 
which allows to identify the Rabi frequency as fi« = 2-f, 
and therefore: 

fi fl = n[j (u) + j 2 (u)} 

= ^ + , (21) 

Now we must solve © using the two-level solutions of 
Eqs.lpll and pOjl. 



. dc n (t) 
1 dt 



fi„ sin {ujt)e Wn0t c a {t)+n nl sin {wt)e iUnlt Cx(t) , 

(22) 



it is easy to show, using (Tl9|) and (|20|) . that 
dc n (t) 



dt 



iOoo/u 



k0 -+- e nk0 



k=—oo 



0* — g* b nfcO* J -)- 



(23) 



where a^ fei = few + cj n j ± a; + 7, = /so; + w ni ± a; - 
7 and (3nki — (fim/4) Jfe(— fijj/a;). We then obtain an 
expression for the coefficients 



c„(t) 



el (fe+l)7r/2 



fc= 


— oo 






a nk0 


b n ko 






a nkl 


Kki ) 



/5 n fco — r b 



4 nfc0 



^nkO 



Pnk\ 



nkl 



■ S n k0 — ibnkl , 



(24) 



where S nkQ = /WVa„ fe0 + VKm ~ l J a nm ~ 1 / h nko) 
and S nkl = Pnki(l/a+ kl + l/b nkl - l/a nkl - l/b+ kl ). 
The average leakage leakage out of the subspace spanned 
by the first two levels is: 
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£ = £Mt)| 2 , 



n=2 



where |c„(i)| 2 means a time average of |c„(i)| 2 . From 
(IM|) we obtain the perturbative result for the average 
leakage as: 



oo oo 



n—2 k— — oo 



1 



1 



( a nfeo) 2 (^nfeo) 2 



« fc l) 2 fccl) 2 

1 1 

+ 



fcl) 2 (d) 2 



^ (5 n A;ocos(fc7r/2) + 5 nfcl siii(fc7r/2) 



— oo / 
/ oo \ 2 

+ ( ^ 5 n fe0 sin (km/ 2) - S nk i cos (kn/2) J . 

\/c=— oo / 

We can simplify the final expression for the leakage 
taking into account that Qa -C u>, since the diagonal 
matrix elements are \Su\ -C 1. In this case the term k = 
of the Bessel functions is dominant in the expansion. We 
obtain 



IV. STRONG DRIVING AND BREAKDOWN 
OF THE TWO-LEVEL APPROXIMATION 

In this section we solve numerically the time-dependent 
Schrodinger equation of the full hamiltonian for the 
Josephson flux qubit, as given by Eq. Qj. As a func- 
tion of time we calcultate the population of the different 
energy eigenstates when a resonant pulse is applied. In 
this way we can evaluate directly how large is the popu- 
lation of the two lowest energy levels and the amount of 
leakage as a function of the pulse strength. 

We first solve numerically the eigenvalue equation of 
Eq. §5$ for fo = 0.5, obtaining the eigenvectors and 
eigenvalues Ei. Then we solve the corresponding time- 
dependent Shcrddinger equation of Eq. Q . 

We start the time evolution at the ground state \*f?o) 
for fo = 0.5. We apply a pulse in the magnetic field, 
/(*) = fo + Sf(t), during a time interval r. After the 
pulse is applied, the wave function has evolved to |^(t)) 
and the qubit returns to / = fo =0.5. We calculate the 
population Pi of the different eigenstates obtained 
at the end of the pulse: Pj(r) = |(*(t)|$,)| 2 . The leak- 
age outside of the quantum computational subspace ex- 
panded by the two lowest eigenstates, \^o), l^i): is then 
obtained as 

oo 

C(j) = l-P (T)-P 1 (T) = Y,Pi(j). 



with 



and 



O 2 O 2 

\Cn(t)\^^ + ^, 



i _ Jg(jWg) ,_ , v , 



Ul n i - U) + 7 UJm + w - 7 



1 



1 



LU ni + LO + 7 LO ni — LO — 7 
1 1 



(W ni - LO + 7) 2 (W ni + LO - 7) 2 
1 1 



2 ' 



(uj ni + LO + 7) 2 (w n i - W — 7) 

Using that Ml n i — 2t: f p aS n i, with S n i — ($„|sin(7r 
4>i — 02 ) I ^o) ; we finally obtain: 



C = 



^ 2 fy^f\s n0 \ 2 , \s nl \ 



h2 \ f 2 

n n=2 V £ «0 



We consider a resonant rf field pulse, 8f(t) = 

(25) fp sin(cj r i) , for < t < t, with the resonant frequency 
huj r = Ei — Eq. In the experimental measurements of 
Rabi oscillations of Ref|6j, pulses of intensity in the range 
5 x 10~ 5 < f p < 5 x 10 -4 were used, as it can be de- 
duced from their data. Here we have calculated the time 
evolution of the Schrodinger equation in a wider range of 
parameters for the pulse strengths: 10~ 5 < f p < 0.1. 

In this section of the work we are going to present re- 
sults for a typical experimental configuration that corre- 
sponds to the work of Chiorescu et alf±. In Fig. 2 we show 
£(t) as a function of the pulse length r for three different 
cases of f p for a — 0.8. We see that £(r) has strong oscil- 
lations as a function of TjQ since high frequencies enter 
into place due to the contribution of several energy levels. 
At f p ~ 10~ 4 , a typical value for experiments, the aver- 
age value of the leakage is very small, C ~ 10~ 7 , showing 
that under a resonant pulse the Josephson flux qubit be- 
haves very closely as a two- level system. In constrast, a 
non-resonant pulse can have a higher leakage for similar 
pulse strengths, as it has been shown in Refllq for a con- 
stant dc pulse. In our results in Fig. 2, we see that for in- 
creasing values of f p the leakage C(t) increases, reaching 
values of 10~ 3 for f p ~ 0.01. A low frequency modulation 
can be clearly observed in Fig. 2(a). This corresponds to 
the Rabi frequency of the approximate two-level system. 
The two-level Rabi frequency increases with f p . Indeed, 

(26) a "Rabi modulation" of the leakage can also be seen in 
Fig. 2(b) at a higher frequency than in Fig. 2(a). On the 
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FIG. 2: Leakage £ as a function of the pulse length t, for an 
rf field pulse of strength f p at the resonant frequency hw r — 
Ei-E . (a) f p = 0.0001; (b) f p = 0.001; (c) f p = 0.01. Time 
is normalized by tj = h/Ej ~ 0.5ps. 



FIG. 4: Populations of the ground state, Po(t) (thick line), 
and the first excited state, Pi(t) (thin line) as a function of 
the pulse length r, for an rf field pulse of strength f p at the 
resonant frequency hw r — E\ — Eq. (a) f p — 0.0001; (b) 
f p = 0.001; (c) f p = 0.01. Time is normalized by tj = H/Ej. 
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FIG. 3: Filled circles: Time averaged leakage £ as a function 
of the pulse strength f p for a — 0.8 for an rf field pulses of 
strength f p at the resonant frequency ftu r = E\ — Eq. Dash- 
dotte line: perturbative approximation of Eq. (|26[) . 



other hand, for f p — 0.01 in Fig. 2(c) the expected Rabi 
frequency is high enough that it can not be distinguished 
from the other high frequencies contributing to the leak- 
age. As we will see below, for f p above this case the 
two-level approximation is no longer adequate. 

To evaluate more quantitatively the effect of strong 
pulses in the amount of leakage, we calculate the time 



averaged leakage £ as a function of the pulse intensity f p . 
We show this result in Fig. 3 for a = 0.8. We observe that 
£ grows quadratically with f p for low pulse strengths. 
At large values of f p the dependence of £ with f p clearly 
departs from this behavior. We compare in Fig. 3 the 
numerical results with the perturbative approximation 
of Eq. ([26|) (summing up to the first 10 levels). We find 
that for f p < /* = 0.02 the perturbative approximation 
is very good. For higher values of the pulse strength, 
the Eq. ([26|) no longer describes the behavior of C(f p ), 
and the the average leakage increases quickly with f p . 
In particular we find that for f p ~ 0.03 the amount of 
leakage is important (i.e., near 10%). From these results, 
we conclude that the two-level approximation can not be 
a good description of the dynamics for pulse strengths 

We show in Fig. 4 the population of the ground state 
Po(t) and the first excited state Pi(t), for the same val- 
ues of f p as in Fig. 2. In Fig. 4(a), for f p — 0.0001, we 
see a clear Rabi oscillation of the populations of the two 
levels of the qubit. A small modulation at high frequen- 
cies is also observed which is due to the perturbation of 
the higher energy levels. In Fig. 4(b), for f p — 0.001, 
the perturbation of the high energy levels is more im- 
portant, but Rabi oscillations can still be distinguished. 
In Fig. 4(c), for f p = 0.01, we see that the behavior of 
Po(t) and P\(t) departs clearly from the simple Rabi 
oscillation scheme and a more complex oscillation with 
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FIG. 5: Triangles: Main oscillation frequency £l oac as a func- 
tion of the pulse strength f p for a = 0.8. Line: Rabi frequency 
Qr as obtained in the rotating-wave approximation, Eq. (|2ip 

two competing frequencies is observed. In this case the 
quantum operation of the qubit, if intended, will be more 
complex since it needs to be based in the knowledge of 
the oscillation pattern. 

We now analize the Fourier power spectrum of the 
populations Po(t) and Pi{t). For example for Pq(t) 
we calculate P(lu) = 1/T\ J Q T exp(W)P 0)^| 2 - We de- 
fine the dominant frequency of the oscillation from the 
maximum of the power spectrum, f2 osc = max tJ P(a;). 
A plot of the obtained f2 osc as a function of the pulse 
strength f p is shown in Fig. 5. We also plot the Rabi fre- 
quency flji given by Eq. (l2ip . which was obtained in the 
rotating-wave approximation. We observe that il osc is 
linear with f p in good agreement with Eq. (j2"Tj) , it i.e. 
f2 osc R3 for f p < = 0.003. For larger values 
of f v , f2 sc departs from the linear dependence of Qr. 
The rotating wave approximation that lead to Eq. ([2~Tj) is 
valid if Qr < u r = (Ei — Eq)/H. Indeed, f2 osc reaches 
this value at fp (we see that for f v — fp — 0.003 we 
have Sl osc /27r » 0.002, while {E x - E )/2n = 0.002107). 
Therefore, for f p > fp we do not expect to find sim- 
ple sinusoidal Rabi oscillations. Instead, more complex 
oscillations are observed, as seen in Fig. 4(c). 

From the above analysis we conclude that for f p > 
f* = 0.003 the possibility of use of the device for quan- 
tum operations becomes more difficult due to the lack of 
Rabi oscillations. This implies that the highest possible 
Rabi frequency that can be obtained is, for a = 0.8, 

f^ Q 727r Q.003aEj/h 5 GHz 

for Ej « (2ttS,)300 GHz, and using tt R 2i:afpEj. In 
the range /* > fp > fp, the two level approximation is 
still a good approximation, since the average leakage is 
relatively small, but its use for quantum operations would 
not be as simple as in the case of f p < fp. For pulse 
strengths f p > f* the intrinsic leakage is very important 



(£ ~ 0.1) and the device can not be treated as a two- level 
system. Moreover, we see in Fig. 5 that at f p ~ /*, the 
dependence of the frequency f2 osc with f p has a drastic 
change, very far apart from a "Rabi regime" . 



V. DEPENDENCE OF THE INTRINSIC 
LEAKAGE WITH CIRCUIT PARAMETERS 

In this section of the paper we are going to study the 
behavior of the leakage as a function of the different pa- 
rameters that define the circuit of the flux qubit. A sim- 
ple argument is that leakage effects should be small if the 
energy difference between the third and the second level 
is much larger than the energy difference between the two 
lowest levels, i.e. if E 2 — E\ ^ E\ — Eq. While this is 
generally correct, the magnitude of the matrix elements 
for the transition rates to higher energy levels can be 
even more important, as we will show in this section. In- 
deed, from the perturbative calculation of Eq. (f2"o) we see 
that, besides an overall factor proportional to fpO 2 , the 
average leakage C depends both on the matrix elements 
S n i, and the factors e n i (which are basically dominated 
by the energy level differences E n — Ei). 

Among the circuit parameters that can be varied either 
by external sources or by circuit design are fo, a and 
Ej/Ec- As we mentioned before the qubit is operated at 
a dc magnetic field near the half-flux quantum, with fo ~ 
1/2. In order to see what happens when we move out of 
the symmetry point fa = 1/2 we plot in Fig. 6(a) the time 
averaged leakage as a function of /o for a fixed value of 
fp. We can see from this figure that the leakage grows 
as we move appart of the symmetry point. This issue 
can be understood with Fig. 6(b). Here we can see that if 
we move out of the symmetry point the distance between 
the two lower levels, Aio = E\ — Eq grows while the third 
level becomes closer to the second one, and thus A21 = 
E-2 — Ei decreases. This basically explains the increase 
of the leakage when moving out of the symmetry point 
fo = 1/2. Further information can be obtained from 
Fig. 6(c) where we show the matrix element IS12I 2 , which 
corresponds to the supercurrent sin(^i — 4>2 — 27r/) taken 
between the first excited level and the second excited 
level. We note that it has a minimum at fo — 0.498. The 
existence of this minimum correlates with the fact that 
the dependence of the leakage with /o has a shoulder near 
this point. 

In principle, the amount of leakage will depend also on 
the a-parameter since this parameter controls the shape 
of the effective potential in the hamiltonian, in particular 
the heigth of the barrier between the two potential min- 
ima. We have calculated numerically the average leak- 
age £ as a function of a for a fixed pulse strengths f p 
at fo = 1/2. This is shown in Fig. 7(a). We see that 
C tends to decrease with increasing a, contrary to what 
the a 2 factor of Eq. (f2"6")) would have suggested. Moreover, 
we see that there is a minimum value for the leakage at 
c*min ~ 0.85. Therefore, there is an optimum value of the 
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FIG. 6: (a) Time averaged leakage £ as a function of fo for 
fixed amplitude f p , Ej/Ec = 35 and a = 0.8. (b) Distance 
between the lower energy levels as a function of fo- (c) Matrix 
element |Si2| 2 as a function of fo. 



FIG. 8: Matrix element 512 as a function of a for different 
values of Ej/Ec- The inset shows the time averaged leak- 
age as a function of Ej/Ec for two different values of the 
amplitude. 
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FIG. 7: (a) Time averaged leakage £ as a function of a for 
fixed amplitude f p , Ej/Ec = 35 and fo = 0.5. (b) Distance 
between the lower energy levels as a function of a. (c) Matrix 
element IS12I 2 as a function of a. 



the third level, the factor that dominates the dependence 
of the leakage with a. 

In Fig. 8 we show the amplitude of the matrix element 
S12 as a function of a for different values of Ej/Eq- We 
observe that the minimum for 1 5*i2 1 2 shown in Fig. 7(c) 
actually corresponds to the fact that 612 crosses zero 
and therefore it vanishes at a particular value of a. Fur- 
thermore, we find that for all the cases of Ej/Ec ana- 
lyzed, the matrix element becomes zero within the range 
of 0.83 < a < 0.89, nearly the same value of a. In the 
inset in Fig. 8 we show how the amount of leakage (for a 
fixed value of f p , a, etc) changes with the ratio Ej/Ec- 
We see that £ does not change significantly, showing a 
slight increase for increasing Ej/Ec- This can be under- 
stood in the sense the dynamics of the system becomes 
more "classical" when increasing Ej/Ec, (the effective 
h decreases), the energy level structure becomes more 
crowded, and thus the effect of higher energy levels will 
tend to be more relevant. 



circuit parameter a for which the leakage will be mimi- 
mum and the two-level approximation more adequate. 

We plot in Fig. 7(b) the difference between the three 
lowest energy levels as a function of a and in Fig. 7(c) 
we plot the matrix element IS12I 2 vs. a. The energy 
distance A 10 strongly decreases with a, while A21 has a 
smooth non-monotonous dependence with a. This later 
result, where A10/A21 decreases with a contributes to 
the general tendency of the leakage to decrease with in- 
creasing a. However, the important result here is that 
the minimum leakage at cv m i n ~ 0.85 is directly corre- 
lated with a minimum in the dependence of IS12I 2 with 
a, as it can be observed in Fig. 7(c). Therefore, it is this 
matrix element, for transitions between the second and 



VI. EFFECTS OF WEAK NOISE ON THE 
INTRINSIC LEAKAGE 

The results of the previous sections have been obtained 
in the ideal case in which the effect of the environment is 
neglected. The aim of this section is to analize how the 
small perturbation of a weak noise can affect the calcu- 
lations of the intrinsic leakage of the previous sections. 
In superconducting qubits various sources of relaxation 
and decoherence are present due to the environment. Re- 
cent experiments in the Josephson flux qubit have shown 
that the dominant source of decoherence is due to 1/f 
noise of the magnetic flux in the SQUID loop4 Here we 
will make the simplification of treating the noise in the 
magnetic flux as a classical noise. In fact, it is usually 
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noise could be in the functional form of the dependence 
of the leakage £(*) for large times (larger than the time 
scale used in the present calculation). 

We can calculate analytically an approximate result 
for the leakage if we assume that the noise intensity A is 
small and that the noise is short-time correlated within 
a small time scale r n , such that A 2 tj <C t„ <C Tr, with 
Tr the period of the Rabi oscillations. In Section III we 
solved the two level system without noise in the RWA 
approximation. Here, we add the effect of noise as a per- 
turbation within the RWA approximation. We now write 
the time dependent perturbation term of the Hamilto- 
nian, W(5f(t)), for small values of f p and small noise 
intensity, as 



FIG. 9: Leakage C as a function of the pulse length r, for an rf 
field pulse of strength f p = 0.00016 at the resonant frequency 
fojJr = Ei — Eq in the presence of white noise of intensity 
A w = 7.10" 6 and 1/f noise with 250 bistable fluctuators and 
intensity Af = 5.10 -5 . Time is normalized by tj = K/Ej ~ 
0.5ps. 



assumed that for time scales much smaller than the en- 
ergy relaxation time, the initial decoherence of the qubit 
can be described with a classical noised Thus, in the 
presence of an rf-drive, we take the magnetic flux as 



/(*) 



1 



f p sm(u)t) + 5f n (t) , 



(27) 



Here, 5f n (t) is the noise source in the magnetic flux, with 
average (Sf n (t)) — and correlations (<5/„(*)<5 /„(*')) — 
A 2 g{t-t'). 

We now solve numerically the time dependent 
Schrodinger equation, Eq.Q, when the rf-pulse is ap- 
plied under the presence of noise, as given by Eg. (1271) . 
We calculate the leakage as in the previous section, but 
now averaging over 50 realizations of the noise. In FigJ§] 
we show the leakage C(t), obtained numerically, as a 
function of the pulse length, for a small value of f p , with 
Gaussian white noise and with 1/f noise. In the first 
case, we have g(t — t') = S(t — *'), and the calculations 
were done for white noise intensity A w ~ 10 _5 . In the 
second case the 1/f noise is defined as the sum of several 
bistable fluctuators as studied, for example, in Ref. HU- 
Here we considered 250 bistable fluctuators with inten- 
sity Af ~ 10~ 5 . As a comparison, in the experiments 
the flux noise intensity is estimated to be of the order 
of A ~ 10 _6 at frequencies of 1Hz4 In the figure we 
can see that, besides the oscillations, there is a general 
linear increase of the leakage as £(*) ~ T^t, which prob- 
ably implies an exponential dependence at long times as 
C(t) s=s 1 — exp(— r^i). As we can see in this figure the 
case of 1// noise and white noise show quite similar be- 
havior in their functional form. (The intensities of the 
two noises shown in Figj9] have been chosen such that 
they result in different leak rates Tl and thus they can 
be distinguished in the plot). The expected difference of 
the realistic case of 1 // noise with a short-time correlated 



W(Sf(t)) « 27rasin(^ + <\> x - fa)(f p sm(wt) + 6f n (t)). 

We then solve equation ([§]) with W n i(t) = 
2iraS n iEj(f p sm(u!t) + 5f n (t)) with the coefficients ob- 
tained from the two level approximation with noise 
(equations (O and ([8])). Integrating equation © and 
keeping only the terms with A 2 and f 2 we obtain 



\Cn(t)\ 2 



= U°) 



(27raA) 5 
h 2 j 



S nO S n iJ (Sl O0 /u])J ({luMP(0) sin 2 (-yt)) 

S n o\ 2 Jo (^Oo/w) 



+ 



7 (P( 7 )+P(- 7 ))t 



+P(0) Bin (2t*)) 

|5m| 2 Jg(nnM 

4 

-P(0) sin(2 7 *) 



7 (P( 7 ) + P(- 7 ))* 



(28) 



where oh (t) is the coefficient obtained in absence of 
noise in equation (|24p . and P(oj) = J^^e 1 ^ 1 g{t)dt is 
the noise spectral density. We see in Eq. (l2"5)) that, be- 
sides the oscillating terms, there are terms that con- 
tribute with an increase of |c„(*)| 2 which is linear in 
f. This will give an overall increase of the leakage as 
£(*) = J2n=2 \ c n(t)\ 2 « lit « 1 - exp(-r L i), which de- 
fines Tl as the leak rate. Typically, in the calculation of 
relaxation and dephasing, a 1/f noise gives a Gaussian 
decay law instcd of an exponential decay. 8 Similarly, one 
could expect that the exponential dependence of the leak- 
age in the presence of short-time correlated noise could 
be modified to a different time dependence in the case of 
1/f noise (possibly a quadratic law), but as we showed 
above, such a difference is not distinguished in our nu- 
merical calculations. 

Then, under the assumption of small short-time corre- 
lated noise, we obtain the leak rate as 



Tr = 



(2naAf 
2h 2 



P 



(29) 
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FIG. 10: Leak rate T L as a function of a for f p — 0.00016 
and white noise strength A = 7 x 10 -6 . Squares: numerical 
result. Circles: perturbative calculation. Inset: Leak rate 
Fl as a function of the white noise strength A for a — 0.8 
and f p — 0.00016. Squares: numerical results. Dashed line: 
perturbative calculation. 
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From the numerical results for the white noise, we can 
make a linear fit of the time dependence of C(t) to obtain 
an estimate of the leak rate T^. In the inset of FigfTUlwe 
show the numerically obtained as a function of the 
white noise intensity A. We compare this result with the 
perturbative calculation of Eq. ([^| . As we can see, the 
agreement is excelent. 

Finally, for a small value of the white noise intensity 
A, we plot the leak rate r^, as a function of the circuit 
parameter a. As we can see in Fig. [TU]the leak rate has a 
minimum near a = 0.85. Therefore we observe that the 
optimum value of a for minimum leakage remains the 
same when a small white noise is added in the system. 
We also compare the numerical result with the perturba- 
tive calculation of Eq. (p?5|) . It shows an overall agreement 
with the existence of a minimum for a = 0.85. We find 
a small systematic difference for a > 0.8. This can be 
attributed to the numerical fit procedure used for esti- 
mating r^, which may need the consideration of larger 
time intervals to improve the reliability of the fit. 



VII. SUMMARY AND DISCUSSION 

We have presented numerical and perturbative calcula- 
tions of the intrinsic leakage in the Josephson flux qubit 
when the system is driven by a rf resonant pulse. As 
a function of the pulse strength f p three regimes have 
been found: (i) For f p < f p , the leakage is very small 
C <C 0.001 and the device shows good Rabi oscillations. 
This is the regime in which the device can be operated 



as a qubit in a simple way. The perturbative calcula- 
tions based on the rotating-wave approximation repro- 
duce very well the numerical results in this case, (ii) For 
fp < f P < fp, the leakage is still small, C < 0.01, and the 
device responds as a two-level system, but the response 
is more complex than simple sinusoidal Rabi oscillations, 
(iii) For f p > /*, the leakage is important, C ~ 0.1 and 
the two-level approximation breaks down. We also find 
that the maximum Rabi frequency that can be achieved 
with this device, for f p = f p , is about £Ir/2tt rj 5GHz. 
This is a factor of 5-10 times larger than what is usually 
achieved in typical experiments, meaning that in princi- 
ple is possible to use stronger pulses and correspondingly 
to further reduce the qubit operation time t op . 

These results refer to the intrinsic leakage, i.e. the 
leakage due to direct transitions to the higher energy 
levels, and therefore the effect of the environment has 
been neglected. The interaction with the environment 
also adds further leakage effects due to indirect transi- 
tions to higher energy levels, see for example Ref ll3lfl4l 
Concerning this point, in Sec. VI we have shown that 
a weak perturbation of the external world, considered as 
a small classical noise, does not change the qualitative 
dependence with circuit parameters of the leakage ob- 
tained in the previous sections for the ideal case. These 
results should be valid as soon as the Rabi period is much 
smaller than the decoherence time, which is the situation 
analyzed in this paper. 

Here we have considered sharp pulses. Several strate- 
gies to reduce leakage varying the pulse shape and the 
type of pulse sequences have been discussed recently*^— 
Applying these strategies the value of the leakage could 
be further reduced from the values obtained here. Ac- 
cording to our results, the use of these strategies will be 
of interest in the regime of f p < f p , when simple Rabi 
oscillations can be observed in the Josephson flux qubit. 

One importan result of this work is the dependence of 
the leakage with the a parameter of the device. We have 
found that there is an optimum value near a ~ 0.85 for 
minimum leakage. This result remains valid when the 
effect of noise of small amplitude in the magnetic flux is 
considered. We have found in this case that the magni- 
tude of the leakage is dominated by the matrix element of 
the supercurrent for transitions between the first and the 
second excited level. When the amplitude of this matrix 
element crosses zero, we find that the leakage is mini- 
mum. Moreover, the Fig. 7 shows that small values of a 
give large leakage, while the range 0.75 < a < 1 is more 
favourable. Of course, the choice of a in the fabrication 
of the device should take into account other factors as 
dependence with a of the gap, the relaxation rate, the 
decoherence time, etc. Interestingly, the experiment of 
Ref|9| shows that the relaxation rate is stronger for small 
values of a, which goes in the same direction as what we 
have observed here for the leakage. 

Recently, Lucero et alS^ have used a procedure (called 
"Ramsey filtering" ) to measure the population of the sec- 
ond excited level in the superconducting phase qubit. It 
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will be interesting if a similar procedure could be imple- 
mented in the Josephson flux qubit. In particular, it will 
be possible to measure experimentally the dependence of 
the leakage with a by using "Ramsey filtering" in the 
recently developed circuit for tuning the gap (i.e. a) of 
the superconducting flux qubit.— ^ 

The results obtained here show quantitatively how the 
multilevel dynamics can become relevant for strong driv- 
ing amplitudes in the flux qubit, going beyond the two- 
level approximation. Indeed, the recent experiments on 
amplitude spectroscopy of Refl27l show how the use of 



strong ac drivings can be turned as a tool to reveal the 
energy level structure of the flux qubit device, now stud- 
ied as a solid-sate artificial atom. 
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